Plot Data
library(ggplot2)
# raw data
ggplot(dataset, aes(x=Bleomycin, y=Counts)) +
theme_bw() +
theme(panel.grid=element_blank(), text = element_text(size=14)) +
geom_smooth(method=lm, formula = y ~ poly(x,2), se=FALSE, aes(colour=siRNA)) +
geom_point(aes(colour=siRNA, shape=Experiment), size=2) +
facet_grid(. ~ genotype) +
xlab(label = "Bleomycin (log10 µM)") +
scale_shape_manual(values=15:20) +
scale_color_manual(values=c("#000000","#FF0000"))

# NormCounts Linear
ggplot(dataset, aes(x=Bleomycin, y=NormCounts, color=siRNA)) +
theme_bw() +
theme(panel.grid=element_blank(), text = element_text(size=14)) +
geom_point(aes(colour=siRNA), size=2) +
geom_smooth(method=lm, formula = y ~ x, se=FALSE) +
facet_grid(. ~ genotype) +
xlab(label = "Bleomycin (log10 µM)") +
scale_color_manual(values=c("#000000","#FF0000"))

# NormCounts2 Linear
ggplot(dataset, aes(x=Bleomycin, y=NormCounts2, color=siRNA)) +
theme_bw() +
theme(panel.grid=element_blank(), text = element_text(size=14)) +
geom_point(aes(colour=siRNA), size=2) +
geom_smooth(method=lm, formula = y ~ x, se=FALSE) +
facet_grid(. ~ genotype) +
xlab(label = "Bleomycin (log10 µM)") +
scale_color_manual(values=c("#000000","#FF0000"))

# NormCounts Quadratic
ggplot(dataset, aes(x=Bleomycin, y=NormCounts, color=siRNA)) +
theme_bw() +
theme(panel.grid=element_blank(), text = element_text(size=14)) +
geom_point(aes(colour=siRNA), size=2) +
geom_smooth(method=lm, formula = y ~ poly(x,2), se=FALSE) +
facet_grid(. ~ genotype) +
xlab(label = "Bleomycin (log10 µM)")+
scale_color_manual(values=c("#000000","#FF0000"))

# NormCounts2 Quadratic
ggplot(dataset, aes(x=Bleomycin, y=NormCounts2, color=siRNA)) +
theme_bw() +
theme(panel.grid=element_blank(), text = element_text(size=14)) +
geom_point(aes(colour=siRNA), size=2) +
geom_smooth(method=lm, formula = y ~ poly(x,2), se=FALSE) +
facet_grid(. ~ genotype) +
xlab(label = "Bleomycin (log10 µM)") +
scale_color_manual(values=c("#000000","#FF0000"))

# NormCounts Cubic
ggplot(dataset, aes(x=Bleomycin, y=NormCounts, color=siRNA)) +
theme_bw() +
theme(panel.grid=element_blank(), text = element_text(size=14)) +
geom_point(aes(colour=siRNA), size=2) +
geom_smooth(method=lm, formula = y ~ poly(x,3), se=FALSE) +
facet_grid(. ~ genotype) +
xlab(label = "Bleomycin (log10 µM)")+
scale_color_manual(values=c("#000000","#FF0000"))

# NormCounts2 Cubic
ggplot(dataset, aes(x=Bleomycin, y=NormCounts2, color=siRNA)) +
theme_bw() +
theme(panel.grid=element_blank(), text = element_text(size=14)) +
geom_point(aes(colour=siRNA), size=2) +
geom_smooth(method=lm, formula = y ~ poly(x,3), se=FALSE) +
facet_grid(. ~ genotype) +
xlab(label = "Bleomycin (log10 µM)") +
scale_color_manual(values=c("#000000","#FF0000"))

library(Cairo)
cairo_pdf("FigureS3B.pdf", width = 5, height = 4, family = "Arial")
ggplot(dataset, aes(x=Bleomycin, y=NormCounts2)) +
theme_bw() +
theme(panel.grid.major=element_blank(), panel.grid.minor=element_blank(),
axis.line = element_line(colour = "black"), text = element_text(size=14),
panel.border = element_blank(), panel.background = element_blank()) +
geom_point(aes(colour = siRNA, shape = genotype), size=1.75) +
geom_smooth(method=lm, formula = y ~ poly(x,2), se=TRUE,
aes(group = GSID,colour = siRNA, linetype = genotype), fill='#DDDDDD', size=0.5) +
#facet_grid(. ~ genotype) +
xlab(label = "Bleomycin (log10 µM)") +
ylab(label = "Normalized Counts") +
scale_color_manual(values=c("#000000","#FF0000")) +
guides(linetype = guide_legend(override.aes= list(color = "#555555")))
dev.off()
## quartz_off_screen
## 2
Models
library(MASS)
library(DHARMa)
library(lme4)
library(lmerTest)
library(bbmle)
Linear formula
fit1 <- lm(Counts ~ Experiment + Bleomycin*siRNA*genotype, data = dataset)
print(summary(fit1))
##
## Call:
## lm(formula = Counts ~ Experiment + Bleomycin * siRNA * genotype,
## data = dataset)
##
## Residuals:
## Min 1Q Median 3Q Max
## -844.27 -280.80 -46.55 227.69 1204.98
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 865.77 253.57 3.414 0.001533 **
## Experimentexp2 660.81 167.46 3.946 0.000331 ***
## Experimentexp3 837.75 167.46 5.003 1.32e-05 ***
## Bleomycin -710.09 224.14 -3.168 0.003025 **
## siRNAsiPARP1 -397.90 331.51 -1.200 0.237463
## genotypeALC1KO 494.49 331.51 1.492 0.144044
## Bleomycin:siRNAsiPARP1 63.74 316.98 0.201 0.841702
## Bleomycin:genotypeALC1KO -518.54 316.98 -1.636 0.110119
## siRNAsiPARP1:genotypeALC1KO -243.48 468.82 -0.519 0.606537
## Bleomycin:siRNAsiPARP1:genotypeALC1KO 376.67 448.27 0.840 0.406014
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 473.6 on 38 degrees of freedom
## Multiple R-squared: 0.7123, Adjusted R-squared: 0.6442
## F-statistic: 10.46 on 9 and 38 DF, p-value: 6.495e-08
cat("AIC: ", AIC(fit1))
## AIC: 738.4095
simres <- simulateResiduals(fittedModel = fit1)
plot(simres)

fit2 <- lm(NormCounts ~ Bleomycin*siRNA*genotype, data = dataset)
print(summary(fit2))
##
## Call:
## lm(formula = NormCounts ~ Bleomycin * siRNA * genotype, data = dataset)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.7392 -0.3245 0.0422 0.3185 0.6835
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.8106 0.2000 9.055 3.12e-11 ***
## Bleomycin -0.9543 0.1912 -4.991 1.22e-05 ***
## siRNAsiPARP1 0.4370 0.2828 1.545 0.1302
## genotypeALC1KO 0.5191 0.2828 1.836 0.0739 .
## Bleomycin:siRNAsiPARP1 -0.5144 0.2704 -1.902 0.0643 .
## Bleomycin:genotypeALC1KO -0.6110 0.2704 -2.260 0.0294 *
## siRNAsiPARP1:genotypeALC1KO -0.5150 0.3999 -1.288 0.2053
## Bleomycin:siRNAsiPARP1:genotypeALC1KO 0.6062 0.3824 1.585 0.1208
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.404 on 40 degrees of freedom
## Multiple R-squared: 0.8402, Adjusted R-squared: 0.8123
## F-statistic: 30.05 on 7 and 40 DF, p-value: 5.166e-14
cat("AIC: ", AIC(fit2))
## AIC: 58.46559
simres <- simulateResiduals(fittedModel = fit2)
plot(simres)

fit3 <- lm(NormCounts2 ~ Bleomycin*siRNA*genotype, data = dataset)
print(summary(fit3))
##
## Call:
## lm(formula = NormCounts2 ~ Bleomycin * siRNA * genotype, data = dataset)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.2685 -0.1196 0.0231 0.1207 0.2482
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.99092 0.07553 13.120 4.47e-16 ***
## Bleomycin -0.52225 0.07222 -7.232 8.94e-09 ***
## siRNAsiPARP1 -0.14330 0.10681 -1.342 0.187
## genotypeALC1KO -0.14489 0.10681 -1.356 0.183
## Bleomycin:siRNAsiPARP1 -0.03162 0.10213 -0.310 0.759
## Bleomycin:genotypeALC1KO -0.04619 0.10213 -0.452 0.654
## siRNAsiPARP1:genotypeALC1KO 0.16411 0.15106 1.086 0.284
## Bleomycin:siRNAsiPARP1:genotypeALC1KO 0.03281 0.14444 0.227 0.821
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1526 on 40 degrees of freedom
## Multiple R-squared: 0.8603, Adjusted R-squared: 0.8359
## F-statistic: 35.19 on 7 and 40 DF, p-value: 3.715e-15
cat("AIC: ", AIC(fit3))
## AIC: -34.99918
simres <- simulateResiduals(fittedModel = fit3)
plot(simres)

fit4 <- lmer(Counts ~ Bleomycin*siRNA*genotype + (1|UID), data = dataset)
print(summary(fit4))
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Counts ~ Bleomycin * siRNA * genotype + (1 | UID)
## Data: dataset
##
## REML criterion at convergence: 634.7
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.7888 -0.4565 -0.1777 0.6047 2.2043
##
## Random effects:
## Groups Name Variance Std.Dev.
## UID (Intercept) 170269 412.6
## Residual 232885 482.6
## Number of obs: 48, groups: UID, 12
##
## Fixed effects:
## Estimate Std. Error df t value
## (Intercept) 1365.29 337.34 16.83 4.047
## Bleomycin -710.09 228.36 32.00 -3.109
## siRNAsiPARP1 -397.90 477.07 16.83 -0.834
## genotypeALC1KO 494.49 477.07 16.83 1.037
## Bleomycin:siRNAsiPARP1 63.74 322.96 32.00 0.197
## Bleomycin:genotypeALC1KO -518.54 322.96 32.00 -1.606
## siRNAsiPARP1:genotypeALC1KO -243.48 674.67 16.83 -0.361
## Bleomycin:siRNAsiPARP1:genotypeALC1KO 376.67 456.73 32.00 0.825
## Pr(>|t|)
## (Intercept) 0.000852 ***
## Bleomycin 0.003920 **
## siRNAsiPARP1 0.415934
## genotypeALC1KO 0.314618
## Bleomycin:siRNAsiPARP1 0.844789
## Bleomycin:genotypeALC1KO 0.118183
## siRNAsiPARP1:genotypeALC1KO 0.722679
## Bleomycin:siRNAsiPARP1:genotypeALC1KO 0.415639
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) Blmycn sRNAsPARP1 gALC1K Bl:RNAPARP1 B:ALC1 sRNAPARP1:
## Bleomycin -0.575
## siRNAsPARP1 -0.707 0.407
## gntypALC1KO -0.707 0.407 0.500
## Bl:RNAPARP1 0.407 -0.707 -0.575 -0.288
## Blmy:ALC1KO 0.407 -0.707 -0.288 -0.575 0.500
## sRNAPARP1:A 0.500 -0.288 -0.707 -0.707 0.407 0.407
## B:RNAPARP1: -0.288 0.500 0.407 0.407 -0.707 -0.707 -0.575
cat("AIC: ", AIC(fit4))
## AIC: 654.7104
simres <- simulateResiduals(fittedModel = fit4)
plot(simres)

Quadratic formula
fit5 <- lm(Counts ~ Experiment + poly(Bleomycin,2)*siRNA*genotype, data = dataset)
print(summary(fit5))
##
## Call:
## lm(formula = Counts ~ Experiment + poly(Bleomycin, 2) * siRNA *
## genotype, data = dataset)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1103.68 -248.93 -34.77 204.19 945.57
##
## Coefficients:
## Estimate Std. Error t value
## (Intercept) 262.6 160.7 1.634
## Experimentexp2 660.8 160.7 4.112
## Experimentexp3 837.8 160.7 5.213
## poly(Bleomycin, 2)1 -3001.1 909.1 -3.301
## poly(Bleomycin, 2)2 180.7 909.1 0.199
## siRNAsiPARP1 -343.8 185.6 -1.852
## genotypeALC1KO 54.0 185.6 0.291
## poly(Bleomycin, 2)1:siRNAsiPARP1 269.4 1285.6 0.210
## poly(Bleomycin, 2)2:siRNAsiPARP1 1012.1 1285.6 0.787
## poly(Bleomycin, 2)1:genotypeALC1KO -2191.6 1285.6 -1.705
## poly(Bleomycin, 2)2:genotypeALC1KO 1616.6 1285.6 1.257
## siRNAsiPARP1:genotypeALC1KO 76.5 262.4 0.292
## poly(Bleomycin, 2)1:siRNAsiPARP1:genotypeALC1KO 1592.0 1818.1 0.876
## poly(Bleomycin, 2)2:siRNAsiPARP1:genotypeALC1KO -1661.6 1818.1 -0.914
## Pr(>|t|)
## (Intercept) 0.111515
## Experimentexp2 0.000235 ***
## Experimentexp3 9.09e-06 ***
## poly(Bleomycin, 2)1 0.002268 **
## poly(Bleomycin, 2)2 0.843612
## siRNAsiPARP1 0.072653 .
## genotypeALC1KO 0.772812
## poly(Bleomycin, 2)1:siRNAsiPARP1 0.835271
## poly(Bleomycin, 2)2:siRNAsiPARP1 0.436588
## poly(Bleomycin, 2)1:genotypeALC1KO 0.097371 .
## poly(Bleomycin, 2)2:genotypeALC1KO 0.217160
## siRNAsiPARP1:genotypeALC1KO 0.772429
## poly(Bleomycin, 2)1:siRNAsiPARP1:genotypeALC1KO 0.387384
## poly(Bleomycin, 2)2:siRNAsiPARP1:genotypeALC1KO 0.367195
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 454.5 on 34 degrees of freedom
## Multiple R-squared: 0.763, Adjusted R-squared: 0.6724
## F-statistic: 8.419 on 13 and 34 DF, p-value: 3.064e-07
cat("AIC: ", AIC(fit5))
## AIC: 737.115
simres <- simulateResiduals(fittedModel = fit5)
plot(simres)

fit6 <- lm(NormCounts ~ poly(Bleomycin,2)*siRNA*genotype, data = dataset)
print(summary(fit6))
##
## Call:
## lm(formula = NormCounts ~ poly(Bleomycin, 2) * siRNA * genotype,
## data = dataset)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.42105 -0.11261 -0.00354 0.08330 0.54949
##
## Coefficients:
## Estimate Std. Error t value
## (Intercept) 1.000e+00 5.461e-02 18.313
## poly(Bleomycin, 2)1 -4.033e+00 3.783e-01 -10.660
## poly(Bleomycin, 2)2 2.759e-01 3.783e-01 0.729
## siRNAsiPARP1 -1.132e-16 7.723e-02 0.000
## genotypeALC1KO -3.754e-16 7.723e-02 0.000
## poly(Bleomycin, 2)1:siRNAsiPARP1 -2.174e+00 5.350e-01 -4.063
## poly(Bleomycin, 2)2:siRNAsiPARP1 2.333e+00 5.350e-01 4.360
## poly(Bleomycin, 2)1:genotypeALC1KO -2.582e+00 5.350e-01 -4.827
## poly(Bleomycin, 2)2:genotypeALC1KO 2.655e+00 5.350e-01 4.961
## siRNAsiPARP1:genotypeALC1KO 2.671e-16 1.092e-01 0.000
## poly(Bleomycin, 2)1:siRNAsiPARP1:genotypeALC1KO 2.562e+00 7.567e-01 3.386
## poly(Bleomycin, 2)2:siRNAsiPARP1:genotypeALC1KO -2.919e+00 7.567e-01 -3.857
## Pr(>|t|)
## (Intercept) < 2e-16 ***
## poly(Bleomycin, 2)1 1.09e-12 ***
## poly(Bleomycin, 2)2 0.470521
## siRNAsiPARP1 1.000000
## genotypeALC1KO 1.000000
## poly(Bleomycin, 2)1:siRNAsiPARP1 0.000251 ***
## poly(Bleomycin, 2)2:siRNAsiPARP1 0.000104 ***
## poly(Bleomycin, 2)1:genotypeALC1KO 2.55e-05 ***
## poly(Bleomycin, 2)2:genotypeALC1KO 1.69e-05 ***
## siRNAsiPARP1:genotypeALC1KO 1.000000
## poly(Bleomycin, 2)1:siRNAsiPARP1:genotypeALC1KO 0.001727 **
## poly(Bleomycin, 2)2:siRNAsiPARP1:genotypeALC1KO 0.000456 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1892 on 36 degrees of freedom
## Multiple R-squared: 0.9685, Adjusted R-squared: 0.9588
## F-statistic: 100.5 on 11 and 36 DF, p-value: < 2.2e-16
cat("AIC: ", AIC(fit6))
## AIC: -11.44408
simres <- simulateResiduals(fittedModel = fit6)
plot(simres)

fit7 <- lm(NormCounts2 ~ poly(Bleomycin,2)*siRNA*genotype, data = dataset)
print(summary(fit7))
##
## Call:
## lm(formula = NormCounts2 ~ poly(Bleomycin, 2) * siRNA * genotype,
## data = dataset)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.152905 -0.045876 -0.001149 0.038348 0.199550
##
## Coefficients:
## Estimate Std. Error t value
## (Intercept) 0.54728 0.02136 25.618
## poly(Bleomycin, 2)1 -2.20725 0.14801 -14.913
## poly(Bleomycin, 2)2 0.15101 0.14801 1.020
## siRNAsiPARP1 -0.17015 0.03021 -5.632
## genotypeALC1KO -0.18412 0.03021 -6.094
## poly(Bleomycin, 2)1:siRNAsiPARP1 -0.13362 0.20931 -0.638
## poly(Bleomycin, 2)2:siRNAsiPARP1 0.83278 0.20931 3.979
## poly(Bleomycin, 2)1:genotypeALC1KO -0.19522 0.20931 -0.933
## poly(Bleomycin, 2)2:genotypeALC1KO 0.91321 0.20931 4.363
## siRNAsiPARP1:genotypeALC1KO 0.19198 0.04273 4.493
## poly(Bleomycin, 2)1:siRNAsiPARP1:genotypeALC1KO 0.13865 0.29601 0.468
## poly(Bleomycin, 2)2:siRNAsiPARP1:genotypeALC1KO -0.99441 0.29601 -3.359
## Pr(>|t|)
## (Intercept) < 2e-16 ***
## poly(Bleomycin, 2)1 < 2e-16 ***
## poly(Bleomycin, 2)2 0.314405
## siRNAsiPARP1 2.16e-06 ***
## genotypeALC1KO 5.20e-07 ***
## poly(Bleomycin, 2)1:siRNAsiPARP1 0.527273
## poly(Bleomycin, 2)2:siRNAsiPARP1 0.000321 ***
## poly(Bleomycin, 2)1:genotypeALC1KO 0.357198
## poly(Bleomycin, 2)2:genotypeALC1KO 0.000103 ***
## siRNAsiPARP1:genotypeALC1KO 6.99e-05 ***
## poly(Bleomycin, 2)1:siRNAsiPARP1:genotypeALC1KO 0.642327
## poly(Bleomycin, 2)2:siRNAsiPARP1:genotypeALC1KO 0.001859 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.074 on 36 degrees of freedom
## Multiple R-squared: 0.9704, Adjusted R-squared: 0.9614
## F-statistic: 107.4 on 11 and 36 DF, p-value: < 2.2e-16
cat("AIC: ", AIC(fit7))
## AIC: -101.5405
simres <- simulateResiduals(fittedModel = fit7)
plot(simres)

fit8 <- lmer(Counts ~ poly(Bleomycin,2)*siRNA*genotype + (1|UID), data = dataset)
print(summary(fit8))
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Counts ~ poly(Bleomycin, 2) * siRNA * genotype + (1 | UID)
## Data: dataset
##
## REML criterion at convergence: 554.3
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.40400 -0.26632 -0.08172 0.47563 1.71356
##
## Random effects:
## Groups Name Variance Std.Dev.
## UID (Intercept) 175352 418.8
## Residual 212556 461.0
## Number of obs: 48, groups: UID, 12
##
## Fixed effects:
## Estimate Std. Error df
## (Intercept) 762.1 276.0 8.0
## poly(Bleomycin, 2)1 -3001.1 922.1 28.0
## poly(Bleomycin, 2)2 180.7 922.1 28.0
## siRNAsiPARP1 -343.8 390.3 8.0
## genotypeALC1KO 54.0 390.3 8.0
## poly(Bleomycin, 2)1:siRNAsiPARP1 269.4 1304.0 28.0
## poly(Bleomycin, 2)2:siRNAsiPARP1 1012.1 1304.0 28.0
## poly(Bleomycin, 2)1:genotypeALC1KO -2191.6 1304.0 28.0
## poly(Bleomycin, 2)2:genotypeALC1KO 1616.6 1304.0 28.0
## siRNAsiPARP1:genotypeALC1KO 76.5 552.0 8.0
## poly(Bleomycin, 2)1:siRNAsiPARP1:genotypeALC1KO 1592.0 1844.2 28.0
## poly(Bleomycin, 2)2:siRNAsiPARP1:genotypeALC1KO -1661.6 1844.2 28.0
## t value Pr(>|t|)
## (Intercept) 2.761 0.02462 *
## poly(Bleomycin, 2)1 -3.255 0.00296 **
## poly(Bleomycin, 2)2 0.196 0.84604
## siRNAsiPARP1 -0.881 0.40413
## genotypeALC1KO 0.138 0.89338
## poly(Bleomycin, 2)1:siRNAsiPARP1 0.207 0.83783
## poly(Bleomycin, 2)2:siRNAsiPARP1 0.776 0.44417
## poly(Bleomycin, 2)1:genotypeALC1KO -1.681 0.10396
## poly(Bleomycin, 2)2:genotypeALC1KO 1.240 0.22538
## siRNAsiPARP1:genotypeALC1KO 0.139 0.89319
## poly(Bleomycin, 2)1:siRNAsiPARP1:genotypeALC1KO 0.863 0.39533
## poly(Bleomycin, 2)2:siRNAsiPARP1:genotypeALC1KO -0.901 0.37526
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) pl(B,2)1 pl(B,2)2 sRNAsPARP1 gALC1K pl(B,2)1:RNAPARP1
## ply(Blm,2)1 0.000
## ply(Blm,2)2 0.000 0.000
## siRNAsPARP1 -0.707 0.000 0.000
## gntypALC1KO -0.707 0.000 0.000 0.500
## pl(B,2)1:RNAPARP1 0.000 -0.707 0.000 0.000 0.000
## pl(B,2)2:RNAPARP1 0.000 0.000 -0.707 0.000 0.000 0.000
## p(B,2)1:ALC 0.000 -0.707 0.000 0.000 0.000 0.500
## p(B,2)2:ALC 0.000 0.000 -0.707 0.000 0.000 0.000
## sRNAPARP1:A 0.500 0.000 0.000 -0.707 -0.707 0.000
## p(B,2)1:RNAPARP1: 0.000 0.500 0.000 0.000 0.000 -0.707
## p(B,2)2:RNAPARP1: 0.000 0.000 0.500 0.000 0.000 0.000
## pl(B,2)2:RNAPARP1 p(B,2)1:A p(B,2)2:A sRNAPARP1:
## ply(Blm,2)1
## ply(Blm,2)2
## siRNAsPARP1
## gntypALC1KO
## pl(B,2)1:RNAPARP1
## pl(B,2)2:RNAPARP1
## p(B,2)1:ALC 0.000
## p(B,2)2:ALC 0.500 0.000
## sRNAPARP1:A 0.000 0.000 0.000
## p(B,2)1:RNAPARP1: 0.000 -0.707 0.000 0.000
## p(B,2)2:RNAPARP1: -0.707 0.000 -0.707 0.000
## p(B,2)1:RNAPARP1:
## ply(Blm,2)1
## ply(Blm,2)2
## siRNAsPARP1
## gntypALC1KO
## pl(B,2)1:RNAPARP1
## pl(B,2)2:RNAPARP1
## p(B,2)1:ALC
## p(B,2)2:ALC
## sRNAPARP1:A
## p(B,2)1:RNAPARP1:
## p(B,2)2:RNAPARP1: 0.000
cat("AIC: ", AIC(fit8))
## AIC: 582.2921
simres <- simulateResiduals(fittedModel = fit8)
plot(simres)

Cubic formula
fit9 <- lm(Counts ~ Experiment + poly(Bleomycin,3)*siRNA*genotype, data = dataset)
print(summary(fit9))
##
## Call:
## lm(formula = Counts ~ Experiment + poly(Bleomycin, 3) * siRNA *
## genotype, data = dataset)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1085.48 -228.62 -45.29 223.66 963.77
##
## Coefficients:
## Estimate Std. Error t value
## (Intercept) 262.56 169.01 1.554
## Experimentexp2 660.81 169.01 3.910
## Experimentexp3 837.75 169.01 4.957
## poly(Bleomycin, 3)1 -3001.14 956.05 -3.139
## poly(Bleomycin, 3)2 180.71 956.05 0.189
## poly(Bleomycin, 3)3 562.81 956.05 0.589
## siRNAsiPARP1 -343.75 195.15 -1.761
## genotypeALC1KO 54.00 195.15 0.277
## poly(Bleomycin, 3)1:siRNAsiPARP1 269.40 1352.06 0.199
## poly(Bleomycin, 3)2:siRNAsiPARP1 1012.10 1352.06 0.749
## poly(Bleomycin, 3)3:siRNAsiPARP1 -865.23 1352.06 -0.640
## poly(Bleomycin, 3)1:genotypeALC1KO -2191.58 1352.06 -1.621
## poly(Bleomycin, 3)2:genotypeALC1KO 1616.58 1352.06 1.196
## poly(Bleomycin, 3)3:genotypeALC1KO -51.65 1352.06 -0.038
## siRNAsiPARP1:genotypeALC1KO 76.50 275.99 0.277
## poly(Bleomycin, 3)1:siRNAsiPARP1:genotypeALC1KO 1591.97 1912.10 0.833
## poly(Bleomycin, 3)2:siRNAsiPARP1:genotypeALC1KO -1661.61 1912.10 -0.869
## poly(Bleomycin, 3)3:siRNAsiPARP1:genotypeALC1KO 272.12 1912.10 0.142
## Pr(>|t|)
## (Intercept) 0.130778
## Experimentexp2 0.000489 ***
## Experimentexp3 2.63e-05 ***
## poly(Bleomycin, 3)1 0.003787 **
## poly(Bleomycin, 3)2 0.851352
## poly(Bleomycin, 3)3 0.560482
## siRNAsiPARP1 0.088354 .
## genotypeALC1KO 0.783903
## poly(Bleomycin, 3)1:siRNAsiPARP1 0.843413
## poly(Bleomycin, 3)2:siRNAsiPARP1 0.459953
## poly(Bleomycin, 3)3:siRNAsiPARP1 0.527071
## poly(Bleomycin, 3)1:genotypeALC1KO 0.115500
## poly(Bleomycin, 3)2:genotypeALC1KO 0.241202
## poly(Bleomycin, 3)3:genotypeALC1KO 0.969779
## siRNAsiPARP1:genotypeALC1KO 0.783538
## poly(Bleomycin, 3)1:siRNAsiPARP1:genotypeALC1KO 0.411663
## poly(Bleomycin, 3)2:siRNAsiPARP1:genotypeALC1KO 0.391747
## poly(Bleomycin, 3)3:siRNAsiPARP1:genotypeALC1KO 0.887783
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 478 on 30 degrees of freedom
## Multiple R-squared: 0.7687, Adjusted R-squared: 0.6376
## F-statistic: 5.864 on 17 and 30 DF, p-value: 1.367e-05
cat("AIC: ", AIC(fit9))
## AIC: 743.9456
simres <- simulateResiduals(fittedModel = fit9)
plot(simres)

fit10 <- lm(NormCounts ~ poly(Bleomycin,3)*siRNA*genotype, data = dataset)
print(summary(fit10))
##
## Call:
## lm(formula = NormCounts ~ poly(Bleomycin, 3) * siRNA * genotype,
## data = dataset)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.42204 -0.04856 -0.01251 0.06747 0.55506
##
## Coefficients:
## Estimate Std. Error t value
## (Intercept) 1.000e+00 5.157e-02 19.391
## poly(Bleomycin, 3)1 -4.033e+00 3.573e-01 -11.288
## poly(Bleomycin, 3)2 2.759e-01 3.573e-01 0.772
## poly(Bleomycin, 3)3 6.525e-01 3.573e-01 1.826
## siRNAsiPARP1 -1.984e-16 7.293e-02 0.000
## genotypeALC1KO -4.092e-16 7.293e-02 0.000
## poly(Bleomycin, 3)1:siRNAsiPARP1 -2.174e+00 5.053e-01 -4.303
## poly(Bleomycin, 3)2:siRNAsiPARP1 2.333e+00 5.053e-01 4.617
## poly(Bleomycin, 3)3:siRNAsiPARP1 -1.425e+00 5.053e-01 -2.821
## poly(Bleomycin, 3)1:genotypeALC1KO -2.582e+00 5.053e-01 -5.111
## poly(Bleomycin, 3)2:genotypeALC1KO 2.655e+00 5.053e-01 5.254
## poly(Bleomycin, 3)3:genotypeALC1KO -6.802e-01 5.053e-01 -1.346
## siRNAsiPARP1:genotypeALC1KO 2.417e-16 1.031e-01 0.000
## poly(Bleomycin, 3)1:siRNAsiPARP1:genotypeALC1KO 2.562e+00 7.146e-01 3.585
## poly(Bleomycin, 3)2:siRNAsiPARP1:genotypeALC1KO -2.919e+00 7.146e-01 -4.084
## poly(Bleomycin, 3)3:siRNAsiPARP1:genotypeALC1KO 1.243e+00 7.146e-01 1.740
## Pr(>|t|)
## (Intercept) < 2e-16 ***
## poly(Bleomycin, 3)1 1.08e-12 ***
## poly(Bleomycin, 3)2 0.445621
## poly(Bleomycin, 3)3 0.077154 .
## siRNAsiPARP1 1.000000
## genotypeALC1KO 1.000000
## poly(Bleomycin, 3)1:siRNAsiPARP1 0.000149 ***
## poly(Bleomycin, 3)2:siRNAsiPARP1 6.03e-05 ***
## poly(Bleomycin, 3)3:siRNAsiPARP1 0.008160 **
## poly(Bleomycin, 3)1:genotypeALC1KO 1.44e-05 ***
## poly(Bleomycin, 3)2:genotypeALC1KO 9.50e-06 ***
## poly(Bleomycin, 3)3:genotypeALC1KO 0.187700
## siRNAsiPARP1:genotypeALC1KO 1.000000
## poly(Bleomycin, 3)1:siRNAsiPARP1:genotypeALC1KO 0.001104 **
## poly(Bleomycin, 3)2:siRNAsiPARP1:genotypeALC1KO 0.000276 ***
## poly(Bleomycin, 3)3:siRNAsiPARP1:genotypeALC1KO 0.091497 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1786 on 32 degrees of freedom
## Multiple R-squared: 0.975, Adjusted R-squared: 0.9633
## F-statistic: 83.23 on 15 and 32 DF, p-value: < 2.2e-16
cat("AIC: ", AIC(fit10))
## AIC: -14.59004
simres <- simulateResiduals(fittedModel = fit10)
plot(simres)

fit11 <- lm(NormCounts2 ~ poly(Bleomycin,3)*siRNA*genotype, data = dataset)
print(summary(fit11))
##
## Call:
## lm(formula = NormCounts2 ~ poly(Bleomycin, 3) * siRNA * genotype,
## data = dataset)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.15326 -0.01990 -0.00463 0.03012 0.20157
##
## Coefficients:
## Estimate Std. Error t value
## (Intercept) 0.54728 0.01926 28.421
## poly(Bleomycin, 3)1 -2.20725 0.13341 -16.545
## poly(Bleomycin, 3)2 0.15101 0.13341 1.132
## poly(Bleomycin, 3)3 0.35710 0.13341 2.677
## siRNAsiPARP1 -0.17015 0.02723 -6.248
## genotypeALC1KO -0.18412 0.02723 -6.761
## poly(Bleomycin, 3)1:siRNAsiPARP1 -0.13362 0.18867 -0.708
## poly(Bleomycin, 3)2:siRNAsiPARP1 0.83278 0.18867 4.414
## poly(Bleomycin, 3)3:siRNAsiPARP1 -0.64854 0.18867 -3.437
## poly(Bleomycin, 3)1:genotypeALC1KO -0.19522 0.18867 -1.035
## poly(Bleomycin, 3)2:genotypeALC1KO 0.91321 0.18867 4.840
## poly(Bleomycin, 3)3:genotypeALC1KO -0.36716 0.18867 -1.946
## siRNAsiPARP1:genotypeALC1KO 0.19198 0.03851 4.985
## poly(Bleomycin, 3)1:siRNAsiPARP1:genotypeALC1KO 0.13865 0.26682 0.520
## poly(Bleomycin, 3)2:siRNAsiPARP1:genotypeALC1KO -0.99441 0.26682 -3.727
## poly(Bleomycin, 3)3:siRNAsiPARP1:genotypeALC1KO 0.57786 0.26682 2.166
## Pr(>|t|)
## (Intercept) < 2e-16 ***
## poly(Bleomycin, 3)1 < 2e-16 ***
## poly(Bleomycin, 3)2 0.266081
## poly(Bleomycin, 3)3 0.011626 *
## siRNAsiPARP1 5.30e-07 ***
## genotypeALC1KO 1.22e-07 ***
## poly(Bleomycin, 3)1:siRNAsiPARP1 0.483934
## poly(Bleomycin, 3)2:siRNAsiPARP1 0.000108 ***
## poly(Bleomycin, 3)3:siRNAsiPARP1 0.001648 **
## poly(Bleomycin, 3)1:genotypeALC1KO 0.308549
## poly(Bleomycin, 3)2:genotypeALC1KO 3.16e-05 ***
## poly(Bleomycin, 3)3:genotypeALC1KO 0.060472 .
## siRNAsiPARP1:genotypeALC1KO 2.08e-05 ***
## poly(Bleomycin, 3)1:siRNAsiPARP1:genotypeALC1KO 0.606887
## poly(Bleomycin, 3)2:siRNAsiPARP1:genotypeALC1KO 0.000749 ***
## poly(Bleomycin, 3)3:siRNAsiPARP1:genotypeALC1KO 0.037891 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.0667 on 32 degrees of freedom
## Multiple R-squared: 0.9787, Adjusted R-squared: 0.9686
## F-statistic: 97.79 on 15 and 32 DF, p-value: < 2.2e-16
cat("AIC: ", AIC(fit11))
## AIC: -109.1627
simres <- simulateResiduals(fittedModel = fit11)
plot(simres)

fit12 <- lmer(Counts ~ poly(Bleomycin,3)*siRNA*genotype + (1|UID), data = dataset)
print(summary(fit12))
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Counts ~ poly(Bleomycin, 3) * siRNA * genotype + (1 | UID)
## Data: dataset
##
## REML criterion at convergence: 491.3
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.26168 -0.23317 0.02957 0.46884 1.68726
##
## Random effects:
## Groups Name Variance Std.Dev.
## UID (Intercept) 168256 410.2
## Residual 240938 490.9
## Number of obs: 48, groups: UID, 12
##
## Fixed effects:
## Estimate Std. Error df
## (Intercept) 762.08 275.98 8.00
## poly(Bleomycin, 3)1 -3001.14 981.71 24.00
## poly(Bleomycin, 3)2 180.71 981.71 24.00
## poly(Bleomycin, 3)3 562.81 981.71 24.00
## siRNAsiPARP1 -343.75 390.29 8.00
## genotypeALC1KO 54.00 390.29 8.00
## poly(Bleomycin, 3)1:siRNAsiPARP1 269.39 1388.35 24.00
## poly(Bleomycin, 3)2:siRNAsiPARP1 1012.10 1388.35 24.00
## poly(Bleomycin, 3)3:siRNAsiPARP1 -865.23 1388.35 24.00
## poly(Bleomycin, 3)1:genotypeALC1KO -2191.58 1388.35 24.00
## poly(Bleomycin, 3)2:genotypeALC1KO 1616.58 1388.35 24.00
## poly(Bleomycin, 3)3:genotypeALC1KO -51.65 1388.35 24.00
## siRNAsiPARP1:genotypeALC1KO 76.50 551.96 8.00
## poly(Bleomycin, 3)1:siRNAsiPARP1:genotypeALC1KO 1591.97 1963.42 24.00
## poly(Bleomycin, 3)2:siRNAsiPARP1:genotypeALC1KO -1661.61 1963.42 24.00
## poly(Bleomycin, 3)3:siRNAsiPARP1:genotypeALC1KO 272.12 1963.42 24.00
## t value Pr(>|t|)
## (Intercept) 2.761 0.02462 *
## poly(Bleomycin, 3)1 -3.057 0.00542 **
## poly(Bleomycin, 3)2 0.184 0.85550
## poly(Bleomycin, 3)3 0.573 0.57178
## siRNAsiPARP1 -0.881 0.40413
## genotypeALC1KO 0.138 0.89338
## poly(Bleomycin, 3)1:siRNAsiPARP1 0.194 0.84778
## poly(Bleomycin, 3)2:siRNAsiPARP1 0.729 0.47306
## poly(Bleomycin, 3)3:siRNAsiPARP1 -0.623 0.53902
## poly(Bleomycin, 3)1:genotypeALC1KO -1.579 0.12753
## poly(Bleomycin, 3)2:genotypeALC1KO 1.164 0.25571
## poly(Bleomycin, 3)3:genotypeALC1KO -0.037 0.97063
## siRNAsiPARP1:genotypeALC1KO 0.139 0.89319
## poly(Bleomycin, 3)1:siRNAsiPARP1:genotypeALC1KO 0.811 0.42544
## poly(Bleomycin, 3)2:siRNAsiPARP1:genotypeALC1KO -0.846 0.40575
## poly(Bleomycin, 3)3:siRNAsiPARP1:genotypeALC1KO 0.139 0.89093
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
cat("AIC: ", AIC(fit12))
## AIC: 527.3356
simres <- simulateResiduals(fittedModel = fit12)
plot(simres)

Final Result
fit <- fit7
output <- coef(summary(fit))
output <- output[grep("Bleomycin", rownames(output)),]
rownames(output) <- gsub("poly\\(|, [1-3]\\)","", rownames(output) )
rownames(output) <- gsub("genotype", paste0(" ",levels(dataset$genotype)[1], " vs. "), rownames(output))
rownames(output)[!(grepl("vs", rownames(output)))] <- paste(rownames(output)[!(grepl("vs", rownames(output)))], levels(dataset$genotype)[1], sep = " in " )
rownames(output) <- gsub("siRNA", paste0(" ",levels(dataset$siRNA)[1], " vs. "), rownames(output))
rownames(output)[!(grepl("vs.*vs| in ", rownames(output)))] <- paste(rownames(output)[!(grepl("vs.*vs| in ", rownames(output)))], levels(dataset$siRNA)[1], sep = " in " )
rownames(output)[!(grepl("vs", rownames(output)))] <- paste(rownames(output)[!(grepl("vs", rownames(output)))], levels(dataset$siRNA)[1], sep = " " )
# suggested result table
kable(output, row.names = T)
| Bleomycin1 in WT siCtrl |
-2.2072495 |
0.1480067 |
-14.9131713 |
0.0000000 |
| Bleomycin2 in WT siCtrl |
0.1510077 |
0.1480067 |
1.0202760 |
0.3144054 |
| Bleomycin1: siCtrl vs. siPARP1 in WT |
-0.1336194 |
0.2093131 |
-0.6383711 |
0.5272726 |
| Bleomycin2: siCtrl vs. siPARP1 in WT |
0.8327816 |
0.2093131 |
3.9786406 |
0.0003209 |
| Bleomycin1: WT vs. ALC1KO in siCtrl |
-0.1952223 |
0.2093131 |
-0.9326808 |
0.3571983 |
| Bleomycin2: WT vs. ALC1KO in siCtrl |
0.9132115 |
0.2093131 |
4.3628969 |
0.0001032 |
| Bleomycin1: siCtrl vs. siPARP1: WT vs. ALC1KO |
0.1386508 |
0.2960134 |
0.4683935 |
0.6423270 |
| Bleomycin2: siCtrl vs. siPARP1: WT vs. ALC1KO |
-0.9944078 |
0.2960134 |
-3.3593333 |
0.0018587 |
write.table(output, file = "FigureS3B_Stats_Ref_WT.txt", quote = F, sep = "\t", row.names = T, col.names = NA)
# re-fit with ALC1KO reference
dataset$genotype <- relevel(dataset$genotype, ref = "ALC1KO")
fit <- lm(NormCounts2 ~ poly(Bleomycin,2)*siRNA*genotype, data = dataset)
output <- coef(summary(fit))
output <- output[grep("Bleomycin", rownames(output)),]
rownames(output) <- gsub("poly\\(|, [1-3]\\)","", rownames(output) )
rownames(output) <- gsub("genotype", paste0(" ",levels(dataset$genotype)[1], " vs. "), rownames(output))
rownames(output)[!(grepl("vs", rownames(output)))] <- paste(rownames(output)[!(grepl("vs", rownames(output)))], levels(dataset$genotype)[1], sep = " in " )
rownames(output) <- gsub("siRNA", paste0(" ",levels(dataset$siRNA)[1], " vs. "), rownames(output))
rownames(output)[!(grepl("vs.*vs| in ", rownames(output)))] <- paste(rownames(output)[!(grepl("vs.*vs| in ", rownames(output)))], levels(dataset$siRNA)[1], sep = " in " )
rownames(output)[!(grepl("vs", rownames(output)))] <- paste(rownames(output)[!(grepl("vs", rownames(output)))], levels(dataset$siRNA)[1], sep = " " )
# suggested result table
kable(output, row.names = T)
| Bleomycin1 in ALC1KO siCtrl |
-2.4024719 |
0.1480067 |
-16.2321812 |
0.0000000 |
| Bleomycin2 in ALC1KO siCtrl |
1.0642192 |
0.1480067 |
7.1903440 |
0.0000000 |
| Bleomycin1: siCtrl vs. siPARP1 in ALC1KO |
0.0050313 |
0.2093131 |
0.0240374 |
0.9809555 |
| Bleomycin2: siCtrl vs. siPARP1 in ALC1KO |
-0.1616262 |
0.2093131 |
-0.7721741 |
0.4450512 |
| Bleomycin1: ALC1KO vs. WT in siCtrl |
0.1952223 |
0.2093131 |
0.9326808 |
0.3571983 |
| Bleomycin2: ALC1KO vs. WT in siCtrl |
-0.9132115 |
0.2093131 |
-4.3628969 |
0.0001032 |
| Bleomycin1: siCtrl vs. siPARP1: ALC1KO vs. WT |
-0.1386508 |
0.2960134 |
-0.4683935 |
0.6423270 |
| Bleomycin2: siCtrl vs. siPARP1: ALC1KO vs. WT |
0.9944078 |
0.2960134 |
3.3593333 |
0.0018587 |
write.table(output, file = "FigureS3B_Stats_Ref_ALC1KO.txt", quote = F, sep = "\t", row.names = T, col.names = NA)